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ABSTRACT 

The formation and growth of supermassive black holes is a key issue to unveil the secrets of galaxy 
formation. In particular, the gravitational recoil produced in the merger of unequal mass black hole 
binaries could have a number of astrophysical implications, such as the ejection of black holes from 
the host galaxy or globular cluster. We present estimates of the recoil velocity that include the effect 
of small eccentricities. The approach is specially suited for the last stage of the merger, where most 
of the emission of linear momentum in gravitational waves takes place. Supplementing our estimates 
with post-Newtonian approximations, we obtain lower and upper bounds that constrain previous recoil 
velocities estimates as well as a best estimate that agrees with numerical simulations in the quasi- 
circular case. For eccentricities e < 0.1, the maximum recoil is found for mass ratios of M1/M2 ^ 0.38 
with a best estimate of ^ 167 (1 -I- e), km s~^ and upper and lower bounds of 79 (1 4- e) km s^^ and 
216 (1 -I- e) km s"'^ respectively. 

Subject headings: black hole physics — cosmology: theory — galaxies: nuclei — gravitation — 
gravitational waves — relativity 



1. INTRODUCTION 

In present hierarchical cold dark matter cosmologies, 
large-scale structures are formed by the merger of small 
(sub-galactic) structures that originated from small pri- 
mordial density perturbations (Press & Schechter 1974). 
In the case of galaxies, there is evidence that most of 
the nearby ones host super-massive black holes (SMBH) 
at their centers (Richstone et al. 1998; Magorrian et al. 
1998) with masses in the range 10^ - lO^Mg. More- 
over, observations have revealed tight relations between 
the SMBH and the bulge of the host galaxy (Ferrarese 
& Merritt 2000), and indicate also that the SMBH mass 
may be determined by the mass of the host dark matter 
halo (Ferrarese 2002). This suggests that there must be 
a deep relation between the formation mechanism of the 
SMBH and the host galaxy which is not yet completely 
understood. It has been suggested (Volonteri et al. 2003) 
that SMBHs grow as a combination of gas accretion and 
mergers with other SMBHs that come together as a re- 
sult of the merger of their host dark matter halos and 
dynamical friction (Begelman et al. 1980). 

The last stages of the merger of SMBHs will be driven 
by the emission of gravitational radiation in the low- 
frequency band. These extremely energetic events will 
be observable by the planned space-based gravitational 
wave antenna LISA (Danzmann 2003; Prince 2003). In 
addition to energy and angular momentum, in the case 
of unequal mass SMBHs, there is a net flux of linear 
momentum carried away form the system by the gravi- 
tational waves (Peres 1962; Bekenstein 1973). Then, mo- 
mentum conservation implies that the final SMBH rem- 
nant will experience a recoil. There may be observational 
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evidence of an ejected SMBH from an ongoing galaxy 
merger (Haehnelt et al. 2006) either from gravitational 
recoil or gravitational slingshot of three o more SMBHs 
in the merger. The knowledge of the magnitude of the 
recoil velocity is crucial to understand the demography 
of SMBHs at the centers of galaxies and in the interstel- 
lar and intergalactic media, and their apparent absence 
in dwarf galaxies and stellar clusters (Madau & Quataert 
2004; Merritt et al. 2004). An estimate of the recoil can 
also be used to constrain theories in which SMBHs grow 
at the center of dark matter halos (Haiman 2004) and to 
estimate SMBH merger rates (Micic et al. 2006). 

In this letter, we present the results of a calculation 
of the recoil velocity based on an approximation scheme 
which is valid when the black holes (BHs) are very close 
to each other, that is, in the last stage of the binary 
merger. It is during this plunge where most of the re- 
coil is accumulated. This contrasts with previous calcu- 
lations that use approximations valid when the BHs are 
well-separated, except for recent fully numerical relativis- 
tic calculations (Baker et al. 2006; Herrmann et al. 2006; 
Gonzalez et al. 2006)). We include the effect of eccentric- 
ity in the recoil velocity since there are some indications 
that SMBH binaries may have not completely circular- 
ized by the time of merger (Aarseth 2003; Armitage & 
Natarajan 2005; Dotti et al. 2006). 

2. RECOIL VELOCITIES AND THE CLOSE LIMIT 
APPROXIMATION 

General relativity predicts that an unequal mass bi- 
nary system in coalescence produces anisotropic emis- 
sion of gravitational waves that carry linear momentum 
away from the binary. As a consequence, the center of 
mass of the system experiences a recoil velocity (Wise- 
man 1992). The first estimates, done by Fitchett (1983) 
using a quasi-Newtonian approach, yielded a maximum 
~ 1480 km s~^ ('"g/^/)^ for a symmetric mass ratio of 
7y = M-^M^/AP ^0.2, where r = 2 GM/c^ is the gravi- 
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tational radius, Mi 2 and M are the individual and total 
masses respectively, and is the orbital separation at 
which gravitational-wave emission ends. Similar results 
up to the innermost stable circular orbit (ISCO) where 
found using a rclativistic perturbative method (Fitch- 
ett & Dctwcilcr 1984). This velocity estimate is largo 
compared to galactic escape velocities. In contrast, the 
first post-Newtonian (PN) calculation (Wiseman 1992) 
for neutron star binaries produced a much lower result: 
t;^ < 1 km . 

Motivated by the astrophysical impact of the gravita- 
tional recoil, there has been recently a number of efforts 
to obtain better estimates of the recoil velocity. Using 
a relativistic perturbative scheme for the description of 
extreme- mass ratio (r; -C 1) binaries and including grav- 
itational back-reaction, Favata et al. (2004) and Merritt 
et al. (2004) obtained maximum recoil at r/ ~ 0.2 with 
167 km s"'^ < Vj, < 622 km s""'^ , depending on the spin 
parameter of the massive black hole (w^ = 456 km s^ for 
the non-spinning case). Recently, Blanchet et al. (2005) 
used second-order PN approximations up to ISCO, fol- 
lowed by an integration of the 2PN linear-momentum 
flux along a plunge geodesic, to obtain a maximum kick 
velocity = 250 ±50 km s^^ for a symmetric mass ratio 
ri = 0.2 . On the other hand, estimates with the effective 
one body approach (Damour & Gopakumar 2006) give 
a maximum recoil in the range 49 — 172 km s^^ for the 
same ij = 0.2. Using the hybrid perturbative-fuUy rel- 
ativistic Lazarus approach, Campanelli (2005) obtained 
also that the maximum recoil takes place at 77 ~ 0.2 
with a magnitude of ^ 250 km s"'^ . Currently, there arc 
very few fully general relativistic numerical simulations 
of kicks. Baker et al. (2006) reported v^. = 105 km s~^ 
for a symmetric mass ratio 77 = 0.24, while Herrmann 
et al. (2006) obtained = 82 km s"^ for rj = 0.23. 
Recently, Gonzalez et al. (2006) report a maximum re- 
coil = 176 ± 11 km s~^ for a symmetric mass ratio of 
T] = 0.195 ± 0.005, where the error bars account only for 
finite differencing Einstein's equations. 

Except for the fully rclativistic approaches, which have 
their own uncertainties (mainly from computational re- 
source limitations), all other approaches become less ac- 
curate the closer the BHs get. Motivated by the fact 
that the main contribution to the recoil comes from the 
last stage of the merger and plunge, as many of the pre- 
viously cited studies suggest, in this letter we compute 
the distribution of kicks using a scheme that works bet- 
ter the closer the orbital separation. This scheme (Price 
& PuUin 1994), known as the close- limit approxima- 
tion (CLA), has been shown to provide accurate results 
for head-on collisions when compared to full numerical 
relativity (Anninos et al. 1995). The CLA is based 
on the idea that the last stages of binary black hole 
(BBH) mergers can be modeled as a single perturbed 
BH. Gravitational recoil was studied in the context of 
the CLA scheme for the case of head-on collisions of 
BHs starting from rest by Andrade & Price (1997). Re- 
cently, we (Sopuerta et al. 2006) extended this analysis 
to the quasi-circular case and found a maximum recoil 
of Vj. ~ 57 km s^^ for rj ~ 0.19 starting from a separa- 
tion of AGM/c^, the maximum separation at which it 
has been estimated that the CLA scheme works (Gleiser 



et al. 1996). This estimate covers the last part of the 
merger, plunge and ringdown. One can complement this 
result with PN estimates from the inspiral phase. If we 
consider the contribution from PN approximations up to 
ISCO, we obtain a lower bound for the maximum re- 
coil of Vj. ~ 80 km s^^ , whereas if we push the PN ap- 
proximations up to the point where we start the CLA 
scheme, we get an upper limit in the maximum recoil 
of ~ 215 km s~^ (it is known that PN and pertur- 
bative schemes overestimate the linear-momentum flux 
in the strong field region). Alternatively, if, instead of 
pushing the PN method towards small separations, we 
push the CLA scheme towards larger separations around 
5 GM/c?, we get a slightly larger upper bound. 

In this letter, we improve and extend these estimates 
for initial configurations corresponding to BBHs in ec- 
centric orbits. The most important point in the imple- 
mentation of the CLA scheme is to establish a corre- 
spondence between the initial configuration representing 
a BBH and a perturbed single BH. This correspondence 
allows us to identify, in the initial BBH configuration, a 
background (a non-spinning Schwarzschild BH) and per- 
turbative multipolar gravitational modes (Sopuerta et al. 
2006). We then evolve these modes by using the machin- 
ery of BH perturbation theory and, from the results of 
the evolution, we can compute the fiuxes of energy, an- 
gular momentum, and linear momentum (-P/j^^) emitted 
by the BBH system during the merger. Then, the recoil 
velocity is given by 

Vr = \\K\\. Vl = -M-' J^' dtP^^. (1) 

where is the time corresponding to the initial sepa- 
ration of the BBH, and tj^ is some time after merger, 

when the gravitational wave emission becomes negligi- 
ble. In our calculations, we use an initial BBH config- 
uration (Brih & Lindquist 1963; Bowen & York 1980) 
that has been shown to accurately represent binaries to 
Newtonian order. They are characterized by four pa- 
rameters (Sopuerta et al. 2006): the total mass of the 
system (A/), the symmetric mass ratio of the binary (ry), 
the initial separation (d), and the initial linear momen- 
tum (P) of each BH. M and r] enter as a scale and input 
parameters, respectively. Given 77, the type of orbit is de- 
termined by the pair (d, P). For a quasi-circular orbit we 
find a relation between P and d by minimizing the grav- 
itational binding energy of the black hole binary. This 
relation turns out to be analogous to the Newtonian one, 
although it has to be interpreted within the framework 
of general relativity. 

In this work, in order to study the effect of eccentricity 
in the gravitational recoil, we take the distance parame- 
ter to be d = a(l — e), where a is a parameter distance 
that plays the role of the semi-major axis and e is the 
eccentricity. Then, the linear momentum parameter is 
given by 

^ = no(a)/^, (2) 

where P^^{a) is the linear momentum parameter for 
quasi- circular configurations (Sopuerta et al. 2006) with 
radius a. 

For our kick estimates, we have kept the semi-major 
axis fixed (a = AGM/c?) and only varied the symmetric 
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mass ratio and eccentricity parameters {f],e). Figure 2 
shows the recoil velocity v^. as a function of the sym- 
metric mass ratio r] for several values of the eccentricity: 
e = — 0.25 . We observe that the maximum recoil oc- 
curs at ry ~ 0.19, consistent with previous results. As we 
increase the eccentricity, the recoil velocity also increases 
by approximately a factor of 1 + e for eccentricities in the 
range e < 0.1. This increase can be qualitatively under- 
stood by expanding equation (2) for small eccentricities, 
which leads to an increase by a factor 1 + e in the linear 
momentum parameter. Thus, it seems that an increase 
in the plunge velocity leads to roughly the same increase 
in the recoil velocity. If we consider the entire range of 
eccentricities that wc considered, the numerical results 
can be fitted to the following two-parameter non-linear 
function 

v^'^^{r), e) = a y/T^ + + P2 V^) 

X (1 + 7i e + 72 + 73 e^) , (3) 

where a (kms~"^), Pi, /Sj, 7i, 72; £^nd 73 are fitting 
parameters given in Table 1 in the CLA row. The 
fit was performed using numerical results in the range 
T] = — 0.25 and e = — 0.5 with an average error 
< 1km s~^ . The fitting function (3) can be shown to 
be functionally equivalent to that presented by Fitchett 
(1983) for small eccentricities, but with different coeffi- 
cients that lead to significantly smaller magnitudes for 
the recoil. 

3. ESTIMATING THE TOTAL RECOIL 

The estimates of the recoil velocities for BBH mergers 
provided by equation (3) and the first row of Table 1 
are incomplete because we have not taken into account 
the contribution from the inspiral (~ 20 km s~^) and 
the beginning of the merger, where fully non-linear com- 
putations are needed. In order to obtain an estimate 
of the total accumulated recoil velocity, we can supple- 
ment our calculations with estimates from approximation 
techniques that are accurate for moderate and large sep- 
arations. In this letter, we use a FN approximation that 
provides the accumulated kick up to some minimum sep- 
aration. For quasi- circular orbits (Blanchet et al. 2005) 
and assuming the angular frequency is given by the New- 
tonian estimate co = y/M/(P, the accumulated recoil up 
to some minimum separations of 4 GM/c^ and 6 GM/ (? 
respectively can be fitted by Eq. (3) with the second and 
third rows of Table 1 respectively. As in Sopuerta et al. 
(2006), upper and lower limits can be obtained by com- 
bining these FN estimates with the CLA calculation. 

Our best estimate consists in adding the FN accumu- 
lated recoil from infinite to a separation of A^GMjc? to 
the CLA estimates, obtained by starting the CLA at that 
separation. However, to compensate for using the FN 
approximation in a regime where it becomes less accu- 
rate, we use a FN expression for the angular velocity 
with errors of 0{v)^ , instead of the Newtonian value used 
in Sopuerta et al. (2006) to obtain an upper limit. The 
results for the quasi-circular case are shown in Fig. 1, 
with error bars coming from the FN approximate error, 
and compared to the most recent full relativistic results 
of Gonzalez et al. (2006) with excellent agreement. 

The main caveat of supplementing the CLA with a 
FN approximation is that presently the latter is only 



available for quasi-circular orbits (Blanchet et al. 2005). 
However, by comparing the contribution from the CLA 
scheme with that provided by the FN approximation for 
the quasi- circular case, we find that the total recoil veloc- 
ity can be accurately described by the CLA contribution 
times a factor: 

^total ^ ^CLA (1 + ^) , (4) 

SO that E can be written in terms of the symmetric mass 
ratio 77 as 

£ = K(1 + Air? + A2r?'), (5) 

where the coefficients k, X^, and A2 depend on the 
minimum distance that we use in the FN approxima- 
tion. When the minimum distance is taken to be 
rf" = 4GM/c2 we obtain: k = 1.357 (0.648), Ai = 
4.418 (5.612), and A2 = 5.33 (23.63) [the estimates in 
parenthesis correspond to using the 2FN angular veloc- 
ity instead of the Keplerian one]; whereas when we take 
it to be rj'" = QGM/c^ the coefficients are: k, = 0.199, 
Al = 4.223, and A2 = 3.926. 

We can then obtain lower and upper bounds as well as a 
best estimate for the total recoil velocities from eccentric 
BBH mergers if we assume that expressions (4) and (5) 
can be extrapolated to the range of eccentricities under 
consideration, e < 0.25. When we use r^"'" = QGM/c^ 
and the Keplerian estimate for the angular velocity, the 
total recoil obtained is a lower limit because we do not 
account for the contribution in the region 4 — QGM/c^. 
On the other hand, if we use rf'" = AGM/(? and the 
Keplerian angular velocity, Eq. (4) provides an upper 
limit since it is well-known that the FN approximation 
overestimate the recoil in the 4 — QGAl/c^ range. Our 
best estimate consists in complementing the CLA with a 
FN estimate using rf = 4 GM/c? and the 2FN angular 
velocity. These estimates, lower and upper bounds and 
best estimate, are shown in Figure 3. 

4. SUMMARY AND DISCUSSION 

We applied the CLA scheme to the computation of re- 
coil velocities for BBH mergers including the effect of 
eccentricity. For small eccentricities, e < 0.1, wc find a 
generic increase in the recoil of the order of 10% with 
respect to the quasi-circular case. This increase is re- 
lated to the fact that, for slightly eccentric orbits, the 
magnitude of the initial velocity of each BH increases 
roughly by a factor of 1 -|- e. Since the CLA scheme is 
valid for separations smaller than the corresponding one 
for the ISCO, we have supplemented our results with es- 
timates obtained with a 2PN approximation, which work 
very well for large separations. Combining appropriately 
both approximation schemes, we have produced lower 
and upper bounds, as well as a best estimate, for the 
recoil velocity. 

Our calculations indicate that the maximum recoil ve- 
locity takes place at a symmetric mass ratio of 77 0.19 
and its magnitude can be as low as ~ 79 — 88 km s^^ 
and as high as ^ 216 — 242 km , with a best esti- 
mate in the range ^ 167 — 187 km s^^ . The variation 
in the upper and lower bounds is due to the eccentricity 
parameter of the orbit. All our estimates can be fit- 
ted by non-linear function of r] and e given by Eq. (3). 
This formula is expected to work well for small eccen- 
tricities, and it reduces to the expression given by Fitch- 
ett (1983), with a significantly smaller magnitude, when 
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TABLE 1 

Non-linear fit parameters for the recoil velocity 



Model 


a 


Pi 


02 


7i 


72 


73 


CLA 


5232 


-2.621 


3.199 


-0.942 


0.808 


-0.405 


PN (rf 


8226 


0.347 


0.083 













1206 


0.1057 


0.05 












^Here rf'" = AGM/c^ and r-j''" = dGM/c^ are the minimum separations up to which the PN approximation is appUed. 



e <C 1. These results also significantly narrow previous 
calculations (Favata et al. 2004; Blanchet et al. 2005; 
Daniour & Gopakumar 2006) and agree remarkably well 
with the most recent full relativistic computations (Gon- 
zalez et al. 2006) in the quasi-circular case. 

One of the main conclusions that can be extracted from 
our and recent estimates of the gravitational recoil from 
BBH mergers is that the main contribution to the recoil 
comes from a narrow interval of separations that include 
the ISCO. Not surprisingly, this coincides with the non- 
linear gravitational regime, where perturbative approxi- 
mation schemes tend to break down. In this sense, the 
attractive feature of our approach is that it uses an ap- 
proximation method, the GLA scheme, valid for small 
separations, and complements the estimate with another 
approximation, the PN scheme, valid for large to moder- 
ate separations. 

As mentioned before. SMBH merger recoil velocities 
are relevant for hierarchical dark matter scenarios for 
structure formation, as well as in the understanding of 
the displacement of nuclear structure in dense galaxies 
and the distortion of X-shaped radio sources. Our results 
do not significantly alter the astrophysical conclusions 
obtained previously by Merritt et al. (2004). Instead, 
our results confirm and sharpen these conclusions, by 



providing a narrower lower and upper bounds for the 
possible recoil velocities. 

Future work can assess the robustness of our estimates 
by studying their dependence on the initial configuration 
as well as the inclusion of higher gravitational multipoles 
in the calculation. Moreover, given the recent advances 
in numerical relativity, it should be possible to compare 
results and to determine with better precision the distri- 
bution of recoil velocities in the parameter space (??, e). 
A very important issue is the effect of spins, which has 
not yet been dealt with enough generality to ascertain 
what is the impact of the additional scales introduced by 
spinning BHs. 
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Recoil velocity in km s as a function of symmetric mass ratio 77 for eccentricities in the range e = — 0.25 . 
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Fig. 3. — Upper and lower limits (dashed and dotted lines) and best estimate (solid lines) to the total recoil velocity in km s~^ as a 
function of symmetric mass ration rj for eccentricities in the range e = — 0.1 . 



